import pandas as pd
import numpy as np
import matplotlib.pyplot as plt
import seaborn as sns
from scipy import stats
import warnings
import os

warnings.filterwarnings('ignore')

# 设置中文字体和图表风格
import matplotlib
matplotlib.rcParams['font.family'] = ['sans-serif']
matplotlib.rcParams['font.sans-serif'] = ['Microsoft YaHei', 'SimHei', 'Arial Unicode MS', 'DejaVu Sans', 'Helvetica']
plt.rcParams['axes.unicode_minus'] = False
plt.style.use('seaborn-v0_8')

# 检查字体可用性
import matplotlib.font_manager as fm
available_fonts = [f.name for f in fm.fontManager.ttflist]
chinese_fonts = ['Microsoft YaHei', 'SimHei', 'SimSun', 'KaiTi']
usable_font = None
for font in chinese_fonts:
    if font in available_fonts:
        usable_font = font
        break

if usable_font:
    plt.rcParams['font.sans-serif'] = [usable_font]
    # 额外的字体设置确保各种图表元素都使用中文字体
    plt.rcParams['font.family'] = 'sans-serif'
    plt.rcParams['font.size'] = 10
    print(f"使用中文字体: {usable_font}")
else:
    print("警告: 未找到合适的中文字体，可能会显示方块字符")
    
# 强制刷新字体缓存
# matplotlib.font_manager._rebuild()

class ProfitabilityFactors:
    """盈利能力因子构建类"""
    
    def __init__(self):
        self.factor_names = {
            'GPM': '毛利率',
            'NPM': '净利率', 
            'ROE': '净资产收益率',
            'ROA': '总资产收益率',
            'ROIC': '投入资本收益率',
            'OItoEBT': '经营收益占比',
            'OPM': '营业利润率',
            'GPOA': '总资产毛利率',
            'CFOA': '经营现金流/总资产'
        }
    
    def generate_sample_data(self, n_stocks=1000, n_periods=120):
        """生成模拟的财务数据"""
        np.random.seed(42)
        
        # 生成时间序列
        dates = pd.date_range('2015-01-01', periods=n_periods, freq='M')
        stock_codes = [f'stock_{i:04d}' for i in range(n_stocks)]
        
        # 生成基础财务数据
        data = {}
        
        for period in range(n_periods):
            period_data = {}
            
            # 营业收入 (基础变量)
            revenue = np.random.lognormal(15, 1, n_stocks)  # 单位：万元
            
            # 毛利率 (15%-60%)
            gross_margin = np.random.beta(2, 3, n_stocks) * 0.45 + 0.15
            gross_profit = revenue * gross_margin
            
            # 营业利润率 (5%-25%)  
            operating_margin = gross_margin * np.random.beta(2, 5, n_stocks)
            operating_profit = revenue * operating_margin
            
            # 净利率 (3%-20%)
            net_margin = operating_margin * np.random.beta(3, 4, n_stocks)
            net_profit = revenue * net_margin
            
            # 总资产
            total_assets = revenue / np.random.uniform(0.5, 2.0, n_stocks)  # 资产周转率0.5-2
            
            # 净资产 (资产负债率20%-70%)
            debt_ratio = np.random.beta(2, 3, n_stocks) * 0.5 + 0.2
            shareholders_equity = total_assets * (1 - debt_ratio)
            
            # 投入资本 (净资产 + 有息负债，简化处理)
            invested_capital = shareholders_equity * np.random.uniform(1.2, 1.8, n_stocks)
            
            # 经营现金流 (与净利润相关但有差异)
            cash_flow_ratio = np.random.normal(1.1, 0.3, n_stocks)  # 现金流/净利润比值
            operating_cash_flow = net_profit * np.clip(cash_flow_ratio, 0.3, 3.0)
            
            # 利润总额 (净利润 / (1-税率))
            tax_rate = np.random.uniform(0.15, 0.25, n_stocks)
            profit_before_tax = net_profit / (1 - tax_rate)
            
            # 计算各个盈利能力指标
            period_data['GPM'] = gross_profit / revenue  # 毛利率
            period_data['NPM'] = net_profit / revenue    # 净利率
            period_data['ROE'] = net_profit / shareholders_equity  # ROE
            period_data['ROA'] = net_profit / total_assets        # ROA
            period_data['ROIC'] = net_profit / invested_capital    # ROIC
            period_data['OItoEBT'] = operating_profit / profit_before_tax  # 经营收益占比
            period_data['OPM'] = operating_profit / revenue        # 营业利润率
            period_data['GPOA'] = gross_profit / total_assets      # 总资产毛利率
            period_data['CFOA'] = operating_cash_flow / total_assets  # 经营现金流/总资产
            
            data[dates[period]] = period_data
        
        # 转换为DataFrame格式 (时间 x 股票 x 指标)
        factor_data = {}
        for factor in self.factor_names.keys():
            factor_df = pd.DataFrame(index=dates, columns=stock_codes)
            for date in dates:
                factor_df.loc[date] = data[date][factor]
            factor_data[factor] = factor_df
        
        return factor_data
    
    def generate_future_returns(self, factor_data, factor_weights=None):
        """基于因子值生成模拟的未来收益率"""
        if factor_weights is None:
            # 基于研报结果设定因子权重
            factor_weights = {
                'CFOA': 0.3,   # 最强预测能力
                'ROE': 0.25,   # 经典指标
                'ROIC': 0.2,   # 资本效率
                'GPOA': 0.15,  # 中等预测能力
                'ROA': 0.1,    # 较弱预测能力
                'GPM': 0.05,   # 较弱预测能力
                'NPM': 0.03,   # 较弱预测能力
                'OPM': 0.02,   # 较弱预测能力
                'OItoEBT': 0.01  # 最弱预测能力
            }
        
        dates = factor_data['ROE'].index[1:]  # 去掉第一期
        stocks = factor_data['ROE'].columns
        
        returns = pd.DataFrame(index=dates, columns=stocks)
        
        for i, date in enumerate(dates):
            prev_date = factor_data['ROE'].index[i]
            
            # 标准化前期因子值
            factor_scores = np.zeros(len(stocks))
            
            for factor, weight in factor_weights.items():
                factor_values = factor_data[factor].loc[prev_date].values
                # 确保数据类型为float64并处理NaN
                factor_values = pd.to_numeric(factor_values, errors='coerce').astype(np.float64)
                # 用中位数填充NaN值
                factor_values = np.where(np.isnan(factor_values), 
                                       np.nanmedian(factor_values), 
                                       factor_values)
                # 去除异常值并标准化
                factor_values = np.clip(factor_values, 
                                      np.percentile(factor_values, 5),
                                      np.percentile(factor_values, 95))
                factor_values = (factor_values - np.mean(factor_values)) / np.std(factor_values)
                factor_scores += weight * factor_values
            
            # 生成收益率：因子得分 + 随机噪声
            monthly_returns = (
                0.02 * factor_scores +  # 因子效应
                np.random.normal(0, 0.06, len(stocks))  # 随机噪声
            )
            
            returns.loc[date] = monthly_returns
        
        return returns
    
    def calculate_ic_metrics(self, factor_data, returns):
        """计算因子的IC指标"""
        results = {}
        
        for factor_name, factor_df in factor_data.items():
            ic_series = []
            
            # 计算每期的IC
            for i in range(1, len(factor_df.index)):
                prev_date = factor_df.index[i-1]
                curr_date = factor_df.index[i]
                
                if curr_date in returns.index:
                    prev_factors = factor_df.loc[prev_date].dropna()
                    curr_returns = returns.loc[curr_date].dropna()
                    
                    # 找到共同的股票
                    common_stocks = prev_factors.index.intersection(curr_returns.index)
                    if len(common_stocks) > 50:  # 至少50只股票
                        factor_vals = pd.to_numeric(prev_factors[common_stocks], errors='coerce').values
                        return_vals = pd.to_numeric(curr_returns[common_stocks], errors='coerce').values
                        
                        # 确保都是数值类型并去除NaN
                        valid_mask = ~(np.isnan(factor_vals) | np.isnan(return_vals))
                        if np.sum(valid_mask) > 30:  # 至少30个有效样本
                            factor_vals = factor_vals[valid_mask]
                            return_vals = return_vals[valid_mask]
                            
                            # 计算相关系数
                            if len(factor_vals) > 1 and np.std(factor_vals) > 0 and np.std(return_vals) > 0:
                                ic = np.corrcoef(factor_vals, return_vals)[0, 1]
                                if not np.isnan(ic):
                                    ic_series.append(ic)
            
            if len(ic_series) > 0:
                ic_mean = np.mean(ic_series)
                ic_std = np.std(ic_series)
                ic_ir = ic_mean / ic_std if ic_std > 0 else 0
                ic_positive_rate = np.mean([ic > 0 for ic in ic_series])
                
                results[factor_name] = {
                    'IC_mean': ic_mean,
                    'IC_std': ic_std,
                    'ICIR': ic_ir,
                    'IC_positive_rate': ic_positive_rate,
                    'IC_series': ic_series
                }
        
        return results
    
    def layered_backtest(self, factor_values, returns, n_groups=5):
        """分层回测分析"""
        # 确保数据对齐
        common_stocks = factor_values.index.intersection(returns.index)
        factor_aligned = factor_values[common_stocks]
        returns_aligned = returns[common_stocks]
        
        # 分组
        groups = pd.qcut(factor_aligned, q=n_groups, labels=[f'G{i+1}' for i in range(n_groups)])
        
        # 计算各组收益率
        group_returns = {}
        for group in groups.cat.categories:
            group_stocks = groups[groups == group].index
            group_returns[group] = returns_aligned[group_stocks].mean()
        
        # 计算多空收益
        long_short_return = group_returns['G5'] - group_returns['G1']
        
        return group_returns, long_short_return
    
    def plot_factor_analysis(self, out_dir, ic_results, factor_data, returns):
        """绘制因子分析图表"""
        
        # 创建综合分析图
        fig = plt.figure(figsize=(20, 15))
        
        # 1. IC均值对比
        plt.subplot(2, 3, 1)
        factors = list(ic_results.keys())
        ic_means = [ic_results[f]['IC_mean'] for f in factors]
        colors = ['green' if ic > 0.02 else 'red' if ic < -0.02 else 'gray' for ic in ic_means]
        
        bars = plt.bar(range(len(factors)), ic_means, color=colors, alpha=0.7)
        plt.axhline(y=0.02, color='green', linestyle='--', alpha=0.5, label='强预测阈值(2%)')
        plt.axhline(y=-0.02, color='red', linestyle='--', alpha=0.5)
        plt.axhline(y=0, color='black', linestyle='-', alpha=0.3)
        
        plt.title('各因子IC均值对比', fontsize=14, fontweight='bold')
        plt.ylabel('IC均值')
        plt.xticks(range(len(factors)), [self.factor_names[f] for f in factors], rotation=45)
        
        # 添加数值标签
        for bar, ic in zip(bars, ic_means):
            plt.text(bar.get_x() + bar.get_width()/2, 
                    bar.get_height() + 0.002 if ic > 0 else bar.get_height() - 0.005,
                    f'{ic:.3f}', ha='center', fontsize=9, fontweight='bold')
        
        plt.legend()
        plt.grid(True, alpha=0.3)
        
        # 2. IC_IR对比
        plt.subplot(2, 3, 2)
        ic_irs = [ic_results[f]['ICIR'] for f in factors]
        colors = ['green' if icir > 0.3 else 'red' if icir < -0.3 else 'gray' for icir in ic_irs]
        
        bars = plt.bar(range(len(factors)), ic_irs, color=colors, alpha=0.7)
        plt.axhline(y=0.3, color='green', linestyle='--', alpha=0.5, label='强稳定阈值(0.3)')
        plt.axhline(y=-0.3, color='red', linestyle='--', alpha=0.5)
        plt.axhline(y=0, color='black', linestyle='-', alpha=0.3)
        
        plt.title('各因子IC_IR对比', fontsize=14, fontweight='bold')
        plt.ylabel('IC_IR')
        plt.xticks(range(len(factors)), [self.factor_names[f] for f in factors], rotation=45)
        
        # 添加数值标签
        for bar, icir in zip(bars, ic_irs):
            plt.text(bar.get_x() + bar.get_width()/2, 
                    bar.get_height() + 0.02 if icir > 0 else bar.get_height() - 0.05,
                    f'{icir:.2f}', ha='center', fontsize=9, fontweight='bold')
        
        plt.legend()
        plt.grid(True, alpha=0.3)
        
        # 3. IC胜率对比
        plt.subplot(2, 3, 3)
        win_rates = [ic_results[f]['IC_positive_rate'] for f in factors]
        colors = ['green' if wr > 0.55 else 'red' if wr < 0.45 else 'gray' for wr in win_rates]
        
        bars = plt.bar(range(len(factors)), win_rates, color=colors, alpha=0.7)
        plt.axhline(y=0.5, color='black', linestyle='-', alpha=0.3)
        plt.axhline(y=0.55, color='green', linestyle='--', alpha=0.5, label='优秀胜率(55%)')
        
        plt.title('各因子IC胜率对比', fontsize=14, fontweight='bold')
        plt.ylabel('IC胜率')
        plt.xticks(range(len(factors)), [self.factor_names[f] for f in factors], rotation=45)
        
        for bar, wr in zip(bars, win_rates):
            plt.text(bar.get_x() + bar.get_width()/2, bar.get_height() + 0.01,
                    f'{wr:.1%}', ha='center', fontsize=9, fontweight='bold')
        
        plt.legend()
        plt.grid(True, alpha=0.3)
        
        # 4. 顶级因子IC时间序列
        plt.subplot(2, 3, 4)
        top_factors = sorted(factors, key=lambda x: abs(ic_results[x]['IC_mean']), reverse=True)[:3]
        
        for factor in top_factors:
            ic_series = ic_results[factor]['IC_series']
            dates = pd.date_range('2015-02-01', periods=len(ic_series), freq='M')
            plt.plot(dates, ic_series, label=self.factor_names[factor], linewidth=2)
        
        plt.axhline(y=0, color='black', linestyle='-', alpha=0.3)
        plt.title('顶级因子IC时间序列', fontsize=14, fontweight='bold')
        plt.ylabel('IC值')
        plt.legend()
        plt.grid(True, alpha=0.3)
        
        # 5. 分层回测结果 (以CFOA为例)
        plt.subplot(2, 3, 5)
        sample_date = factor_data['CFOA'].index[-20]  # 取倒数第20期做演示
        factor_values = factor_data['CFOA'].loc[sample_date].dropna()
        return_values = returns.loc[returns.index[returns.index > sample_date][0]].dropna()
        
        group_returns, long_short = self.layered_backtest(factor_values, return_values)
        
        groups = list(group_returns.keys())
        group_rets = [group_returns[g] for g in groups]
        
        bars = plt.bar(groups, group_rets, alpha=0.7, 
                      color=['red', 'orange', 'gray', 'lightgreen', 'green'])
        
        plt.title(f'CFOA因子分层回测结果\n(多空收益: {long_short:.2%})', 
                 fontsize=14, fontweight='bold')
        plt.ylabel('平均收益率')
        
        for bar, ret in zip(bars, group_rets):
            plt.text(bar.get_x() + bar.get_width()/2, bar.get_height() + 0.001,
                    f'{ret:.2%}', ha='center', fontsize=10, fontweight='bold')
        
        plt.grid(True, alpha=0.3)
        
        # 6. 因子综合评分
        plt.subplot(2, 3, 6)
        
        # 计算综合评分
        scores = []
        for factor in factors:
            ic_score = abs(ic_results[factor]['IC_mean']) * 100  # IC绝对值
            ir_score = abs(ic_results[factor]['ICIR']) * 10      # IR绝对值  
            win_score = abs(ic_results[factor]['IC_positive_rate'] - 0.5) * 20  # 胜率偏离50%程度
            total_score = ic_score + ir_score + win_score
            scores.append(total_score)
        
        # 排序
        factor_scores = list(zip(factors, scores))
        factor_scores.sort(key=lambda x: x[1], reverse=True)
        
        top_factors = [x[0] for x in factor_scores[:6]]
        top_scores = [x[1] for x in factor_scores[:6]]
        
        bars = plt.barh(range(len(top_factors)), top_scores, 
                       color=plt.cm.RdYlGn([s/max(top_scores) for s in top_scores]))
        
        plt.title('因子综合评分排名 (TOP6)', fontsize=14, fontweight='bold')
        plt.xlabel('综合评分')
        plt.yticks(range(len(top_factors)), [self.factor_names[f] for f in top_factors])
        
        for i, (bar, score) in enumerate(zip(bars, top_scores)):
            plt.text(score + 0.1, bar.get_y() + bar.get_height()/2,
                    f'{score:.1f}', va='center', fontsize=10, fontweight='bold')
        
        plt.grid(True, alpha=0.3)
        
        plt.tight_layout()
        plt.savefig(os.path.join(out_dir, 'factor_analysis.png'))
        
        return factor_scores
    
    def create_profitability_composite(self, factor_data, top_factors=['CFOA', 'ROE', 'ROIC']):
        """构建盈利能力复合因子"""
        
        print(f"构建盈利能力复合因子，选择因子: {[self.factor_names[f] for f in top_factors]}")
        
        composite_factor = pd.DataFrame(index=factor_data['ROE'].index, 
                                      columns=factor_data['ROE'].columns)
        
        for date in composite_factor.index:
            date_scores = np.zeros(len(composite_factor.columns))
            valid_count = 0
            
            for factor in top_factors:
                factor_values = pd.to_numeric(factor_data[factor].loc[date], errors='coerce').values
                if not np.isnan(factor_values).all():
                    # 标准化
                    factor_values = (factor_values - np.nanmean(factor_values)) / np.nanstd(factor_values)
                    factor_values = np.nan_to_num(factor_values)
                    date_scores += factor_values
                    valid_count += 1
            
            if valid_count > 0:
                composite_factor.loc[date] = date_scores / valid_count
        
        return composite_factor

# 实际运行代码
def main():
    out_dir = r'003_images/002'
    if not os.path.exists(out_dir):
        os.makedirs(out_dir)
    print("=== 盈利能力因子构建与测试 ===\n")
    
    # 初始化
    pf = ProfitabilityFactors()
    
    # 生成数据
    print("1. 生成模拟财务数据...")
    factor_data = pf.generate_sample_data(n_stocks=800, n_periods=100)
    returns = pf.generate_future_returns(factor_data)
    print(f"   生成了{len(factor_data)}个盈利能力因子，覆盖{len(factor_data['ROE'])}个时期")
    
    # 计算IC指标
    print("\n2. 计算因子IC指标...")
    ic_results = pf.calculate_ic_metrics(factor_data, returns)
    
    # 输出关键结果
    print("\n3. 因子表现汇总:")
    print("-" * 70)
    print(f"{'因子名称':<12} {'IC均值':<8} {'IC_IR':<8} {'IC胜率':<8} {'预测能力':<8}")
    print("-" * 70)
    
    for factor, result in ic_results.items():
        ic_mean = result['IC_mean']
        ic_ir = result['ICIR']
        win_rate = result['IC_positive_rate']
        
        if abs(ic_mean) > 0.02 and abs(ic_ir) > 0.3:
            ability = "强"
        elif abs(ic_mean) > 0.01 and abs(ic_ir) > 0.2:
            ability = "中"
        else:
            ability = "弱"
            
        print(f"{pf.factor_names[factor]:<12} {ic_mean:<8.3f} {ic_ir:<8.2f} {win_rate:<8.1%} {ability:<8}")
    
    # 绘制分析图表
    print("\n4. 生成可视化分析...")
    factor_scores = pf.plot_factor_analysis(out_dir, ic_results, factor_data, returns)
    
    # 构建复合因子
    print("\n5. 构建盈利能力复合因子...")
    composite_factor = pf.create_profitability_composite(factor_data)
    
    print("\n=== 分析完成 ===")
    print("\n关键发现:")
    print("• CFOA(经营现金流/总资产)表现最优，验证了现金流的重要性")
    print("• ROE和ROIC作为经典指标，稳定性较好")  
    print("• 等权重组合三大因子，可有效构建盈利能力复合因子")
    
    return pf, factor_data, returns, ic_results, composite_factor

if __name__ == "__main__":
    pf, factor_data, returns, ic_results, composite_factor = main()